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1.  INTRODUCTION 


In  an  analysis  of  seismic  data  made  for  the  purpose 
of  discriminating  between  earthquakes  and  explosions,  a 
common  practice  (Booker  and  Mitronovas  (1964),  Capon, 
et.  al.  (1967),  Ericsson  (1970),  and  Shumway  and  Blandford 
(1970))  is  to  use  the  mean  of  a  seismic  discriminant 
measured  at  several  stations  as  the  basic  observation. 

For  example,  surface  wave  and  body  wave  pairs  plotted  in 
two  dimensions  may  actually  be  mean  values  deduced  from 
differing  numbers  or  configurations  of  stations.  The 
discriminant  analysis  which  results  ignores  the  fact  that 
the  measurements  may  have  different  means  due  to  differ¬ 
ing  station  configurations  and  will  surely  have  different 
variances  because  of  the  different  sample  sizes  over 
which  each  mean  has  been  determined.  For  these  reasons, 
one  would  prefer  to  work  with  the  values  observed  at 
each  of  the  stations  which  record  the  event. 

This  report  will  describe  and  evaluate  a  method  for 
handling  the  discriminant  analysis  when  observations  are 
made  at  more  than  one  station.  A  linear  criterion  which 
is  an  extension  of  the  standard  discriminant  function 
will  be  considered.  Evaluation  will  be  in  terms 
of  standard  operating  characteristic  curves  (termed 
"identification  curves"  by  Ericsson,  (1970)  which 
measure  the  predicted  explosion  detection  capability 
as  a  function  of  the  false  alarm  rate.  Ericsson,  (1970) 
contains  a  clearly  written  explanation  of  the  political 
ramifications  leading  to  the  use  of  identification  curves  as 
well  as  an  excellent  review  of  the  seismic  discrimination 
literature. 
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2.  LINEAR  DISCRIMINANT  ANALYSIS 
WITH  MULTIPLE  OBSERVATIONS 


In  the  standard  situation  a  p  x  1  vector  of  proposed 
discriminants,  X*  =  (X^,...,X  ),  is  observed  and  we  wish 
to  determine  whether  the  vector  came  from  a  population 
of  earthquakes  (hypothesis  H^)  or  a  population  of  explo¬ 
sions  (hypothesis  H2).  In  general,  multivariable 
normality  of  the  vector  X  is  assumed  with  only  the  mean 
value  vectors  for  the  two  populations  allowed  to  vary. 

The  effect  of  variations  likely  to  be  observed  in  the 
covariance  matrices  seems  to  be  small  for  seismic  data 
(Shumway  and  Blandford,  1970).  In  the  case  where  Ng 
stations  record  observations  from  the  same  event,  we 
would  like  to  assign  each  station  a  different  mean  value. 
Then,  one  has  available  a  sample  of  Ng  vector  discriminants 

(X. ,  j  =  1,  2,...,N  )  for  classifying  the  events  as 
J  5 

belonging  to  either  the  earthquake  or  explosion  group. 

Under  (earthquakes)  ,  (X.. ,  j  =  1,  2,...Ns)  are  Ns 
independent  normal  vectors  with  mean  (V^ j ,  j  =  1,...,NS) 
and  covariance  matrix  l  while  under  H2  they  are  inde¬ 
pendent  normal  vectors  with  means  (V~.,  j  =  1,...,N  ) 

^  J  s 

and  identical  covariance  matrix  £.  An  application  of  the 
Neyman-Pearson  argument  to  this  case  yields  the  dis¬ 
criminant  function: 

Ns  N 

^  *  •  •  •  >  )  1  ^ .  x .  -  £  (  ct  ) . 

S  j=l  J  J  j=l  J 


where 


and 


(CT).  -  1/2^.  +  V2j)«  I'1  CY.ij  "  — 2 j ^  ^ 

Equations  (1)  through  (3)  demonstrate  that  the  discrimi¬ 
nant  function  for  observations  from  more  than  one  station 
is  just  the  sum  of  the  ordinary  discriminant  functions 
evaluated  at  the  station  means.  It  may  be  that  general 
linear  procedures  of  the  form 

N 

a(X^,...,X^  )  £  £*•  X. 

Ns  j-1  3  J  (4) 

might  be  useful  discriminants.  An  example  of  a 

linear  function  which  is  not  necessarily  the  linear 

discriminant  function  would  be  M  .  -  1.9  m,  .  where  the 

s  j  b  j 

vector  Xj  =  (Msj ,  m^j)  is  a  surface  wave  body  wave 
pair  and  aj  =  (1,  -1.9).  We  note  parenthetically  th.t 
for  Ns  =  1,  or  for  equal  station  means,  equations  (1) 
through  (3)  reduce  to  the  usual  equations  for  the 
linear  discriminant  function,  so  that  the  procedure 
is  equivalent  to  discrimination  using  network  averages. 

In  order  to  apply  equations  (1)  through  (3)  we 
must  either  know  Yx . ,  ,  J,  j  =  l,...,Ns,  or  obtain 

consistent  estimators  for  these  terms. 

In  Section  4  we  discuss  the  usual  multivariate 
theory  as  it  applies  to  estimating  the  generalized 
means  and  covariance  matrix. 
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3.  THE  IDENTIFICATION  CURVE  AS  A 
MEASURE  OF  PERFORMANCE 

For  discriminant  functions  of  the  form  (1)  an 
expected  future  performance  can  be  calculated  using  the 
result  that  d(X^,...,X^  )  is  a  normal  random  variable 

with  mean  £  under  (earthquake) ,  and  - jtD?  under  H 2 

(explosion),  and  a  variance  2£  D.  under  both  hypotheses, 
with  J  3 


Dj  -  1/2 (VXj 


^ij  ■  v2j) 


(5) 


In  general  we  wish  to  control  the  probability  of  a  false 
alarm.  A  false  alarm  occurs  when  an  earthquake  is 
incorrectly  called  an  explosion.  That  is  when  Hj  is  true 
and  the  value  of  the  discriminant  function  dfXj,...,^  ) 
is  low  enough  to  cause  the  observation  to  be  classified 
with  the  explosion  population.  If  a  is  the  false  alarm 
probability  and  C(a)  is  a  constant  depending  only  on  a 
we  accept  (explosion)  if 


»  •••»  ^  C  (a)  . 

s 

For  a  given  false  alarm  probability 
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a  ~  ^rH-^  ^  i »  •  •  •  i  C(oi)) 


=  $ 


(C(a)  "  ^Di 

2lB] 


(6) 


with  $  denoting  the  cumulative  normal  distribution 
function,  equation  (6)  may  be  solved  for  C(a)  yielding 


c(°0  =  (I  D?)  (1  +  2  $"1(a)) 

j  J 


(7) 


_2 

where  $  denotes  the  inverse  of  the  cumulative  normal 
distribution  function.  The  probability  of  correctly 
identifying  the  explosion. 


Pd(°0  PrH^  (dQLi»  •••  >£Ns)  1c(°0) 


=  $ 


C 

(- 


(a) 


(8) 


can  be  calculated  as  a  function  of  the  false  alarm  rate  ct 
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We  note  that  for  a  general  linear  function  of  the  form 

2  2 

(4)  with  means  y  ^ ,  ^2  <  and  variances  o^,  a^, 
equations  (7)  and  (8)  become 

C (a)  =  yx  +  cx  $“1(a)  (93 


and 


C(a)  -  y y 
pd(«)  -  *  C— a- - - 


(10) 


If  y£  >  y^  the  equations  become 
C '  (a)  =  yx  +  ax  «f_1  (1-a) 


(11) 


and 


C(a)  -  y- 

p.(a)  -  l  -  *  ( - — -£) 


(12) 


A  plot  of  P^(a),  Pj(a)  against  a  gives  a  measure  of 
explosion  detection  capability  and  has  been  termed  an 
identification  curve  by  Ericsson  (1970). 

A  hypothetical  identification  curve  is  shown  in 
Figure  1.  In  general,  Ericsson  postulates  that  an  accept¬ 
able  operating  region  is  determined  by  keeping  the  false 
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DETECTION  PhJBABIUTY  Pd  (8) 


Figure  1.  A  hypothetical  identification  curve  showing  the  region 
which  combines  an  acceptable  detection  probability  Pd>.10  with 
an  acceptable  risk  a  <  10"3. 


alarm  probability  at  1  x  10  ^  or  less  while  maintaining 
the  desired  detection  probability  at  .10  or  greater. 

Figure  1  gives  a  geometric  interpretation  to  the  politically 
acceptable  region  along  with  a  hypothetical  identification 
curve  showing  that  a  maximum  of  35$  of  the  explosions 

could  be  detected  within  the  acceptable  false  alarm 

_  3 

probability  of  1  x  10 

Figure  1  suggests  that  two  proposed  procedures 
leading  to  normally  distributed  discriminants  can  be 
compared  by  examining  their  respective  identification 
curves.  A  procedure  is  superior  to  another  if  its  detec¬ 
tion  curve  is  located  to  the  right  of  the  other  procedure 
over  the  region  of  interest. 


-7- 


4.  ESTIMATION  OF  Till:  STATION  MEANS  AMD  COVARIANCE 
MATRIX  FOR  EARTHQUAKES  AND  EXPLOSIONS 


In  order  to  estimate  the  unknown  parameters  j , 

-Zj ’  j  =  1 . Ns  we  suppose  that  multiple  station 

observations  are  available  for  N^  earthquakes  and  N£ 
explosions  on  possible  discriminants .  Denote  the  k'th 
p  x  1  vector  of  discriminants  measured  at  station  j  for 
an  event  of  type  i  by  so  that  under  (earthquake) 

has  a  normal  distribution  with  mean  vector  and 
covariance  matrix  £  while  under  (explosion)  £2jk  ^as 
a  mean  vector  V 2 j  and  covariance  matrix  l,  j  =  1,  2,...,.); 
k  =  1,  2,...,  Njj  is  the  number  of  observations  ' 

measured  by  the  j'th  sta:ion  for  an  earthquake  and  N2j 
is  the  number  of  observations  made  by  the  j'th  station 
for  an  explosion. 

In  this  case  it  is  well  known  that  the  maximum 
likelihood  estimators  for  the  mean  vectors  are 


V.  .  = 

ij 


NU  Jx  ’  1  =  1 


=  1»  2,  j  =  1 . J 


and  that  an  unbiased  estimate  for  the  covariance  matrix 
which  is  proportional  to  the  maximum  likelihood  estimator 


*  ■  ijk  ’  ’  V' 
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In  practice  the  computations  are  rather  involved  for 
the  above  expression  and  one  may  wish  to  try  different 
models  and  test  various  hypotheses  about  the  vectors  Xjj* 
Therefore,  it  is  convenient  to  formulate  all  models  in 
terms  of  the  multivariate  linear  hypothesis  structure 
described  below. 

Let  ,  k  =  1,...,N  be  p  x  1  independently  dis¬ 
tributed  normal  vectors  with  means  B'xj.  and  covariance 
matrix  Here  R  is  a  q  x  p  matrix  of  unknown  parameters, 
and  x^,  k  =  1,...,N  are  fixed  known  q  x  1  vectors. 

Then  the  Nxp  matrix  Y  =  (y.i»»..,y  )'  has  nean  value 

N 

XB,  where  X  *  (x^  , . . .  ,x^)  ' .  For  this  ca  e  the  maximum 
likelihood  estimace  for  B  is 


B  =  (X'X)'1  X ' Y 


(15) 


while  an  unbiased  estimate  proportional  to  the  maximum 
likelihood  estimator  for  ][  is 


1  =  (N-q)’1  (Y-XB)'(Y-XB) 


(16) 


For  this  setup,  the  program  from  Dixon  (1970)  yields  B  and 
E,  as  well  as  the  results  of  testing  hypotheses  of  the  form 

ABC*  *  D  (17) 

where  A(r  x  q) ,  B(q  x  s),  and  D(r  x  s)  are  matricc  specified 
to  generate  a  designated  hypothesis.  The  matrix  C  must  be 
(s  x  s) .  We  illustrate  the  procedure  in  the  next  section. 
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5.  liXAMI’Ll: 


I  he  procedure  implied  by  the  theoretical  approach 
in  Sections  2,  3,  and  4  will  be  illustrated  using  20 
earthquakes  and  9  explosions  observed  at  various  sub • 
combinations  of  six  recording  stations.  The  observations 
shown  in  Fable  I  arc  on  body  wave  (m^)  and  surface  wave 
(Mg)  magnitudes.  We  note  that  most  events  arc  not 
recorded  at  all  stations.  Hence,  it  would  be  useful  to 
estimate  an  expected  explosion  detection  probability  for 
each  sub-combination  of  stations  which  can  be  expected 
to  appear  as  an  observation. 

We  assume  as  in  the  previous  section  that  the  k'th 
2x1  vector  =  ^Msijk»Mbijk^  '  pair  at  sta*ion 

j  =  1,  2,... ,6  for  event  i  =  1,  2  (earthquake  or  explo¬ 
sion)  has  the  2x1  mean  vector: 


(18) 


where  components  denote  respectively  the  theoretical 
mean  of  the  Mg  and  nb  measurements  at  station  j  for  an 
event  of  type  i.  In  general,  the  estimates  for  V.  .  and  T 

— 1  J  L 

arc  most  important  as  they  arc  needed  for  the  multiple 
station  discriminant  function  given  in  Section  2. 
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TABU:  I 


I’alrs  for  29  Lvents 
Recrrded  at  6  LRSM  Stations 


Stations 


Larthqunkcs  (1) 

■b 

NP 

M 

s 

mb 

WII 

Ms 

% 

PC 

f,s 

mb 

KT 

Ms 

"b 

IIN 

M 

5 

MB 

KN 

Ms 

LQI 

6.25 

4.30 

5.07 

4.20 

5.99 

4.94 

LQ2 

5.77 

4.23 

4  80 

3.81 

4.86 

3.75 

LQ3 

6.26 

6.33 

6.33 

6.39 

6.35 

6.19 

LQ4 

6.24 

4.  15 

S.  24 

4.70 

4.61 

4.62 

hqs 

6.46 

6.17 

6.40 

6.29 

6.04 

6.57 

S.62 

6.46 

5.29 

6.46 

LQ6 

6.21 

4.59 

4.86 

4.52 

5.01 

4.36 

4.98 

4.20 

LQ7 

5.45 

5.01 

5.07 

S.  38 

4.98 

4.91 

LQ8 

5.56 

4.41 

4.52 

4.74 

4.61 

4.28 

4.47 

4.56 

4.61 

4.51 

LQ9 

5.89 

5.36 

S.  22 

5.67 

S.  30 

5.42 

5.4  3 

5.48 

S .  36 

5.62 

4.96 

5.39 

LQ10 

5.44 

4.93 

4.66 

4.30 

4.94 

4.21 

LQI  1 

5. 22 

4.40 

5.30 

4.22 

5.23 

4.62 

5.31 

4.40 

i:qi2 

5.85 

4.69 

5.18 

4.73 

5.40 

4.50 

5.31 

4.73 

5.24 

4.77 

4.62 

4.80 

i:gu 

5.37 

3.58 

4.62 

3.67 

5.13 

3.58 

4.78 

3  81 

EQI4 

5.81 

3.86 

4.92 

3.90 

4.70 

4.08 

5.18 

4.10 

4.95 

4.01 

LQI  5 

6.05 

4.51 

4.98 

4.93 

4.85 

4.84 

4.89 

4.?S 

4.80 

4.52 

4.78 

4.75 

LQ  1 6 

5.90 

4.86 

5.06 

5.05 

4.94 

4.90 

4.83 

4.81 

4.90 

4.84 

4.85 

4.95 

LQI  7 

5.81 

4.94 

4.95 

4.75 

5.15 

4.95 

4.72 

4.57 

LQI  8 

4.51 

4.36 

4.54 

3.95 

4.63 

3.92 

LQI  9 

5.98 

5.00 

4.83 

5.17 

4.72 

4.91 

4.67 

5.02 

5.21 

4.89 

4.60 

4.65 

LQ20 

5.51 

5.40 

5.48 

5.35 

5.23 

5.  15 

5.78 

5.33 

5.17 

5.11 

I’rosuraod  ... 

Explosions 

LX1 

6. 19 

4.10 

5.81 

4.58 

6.12 

4.31 

LX  2 

6.13 

4.33 

6.08 

4.45 

5.99 

4.40 

5.75 

4.25 

5.93 

4.99 

LX  3 

6.25 

4.32 

6.08 

4.45 

5.99 

4.40 

5.73 

4.26 

LX  4 

6.02 

3.R4 

6.25 

3.76 

5.84 

4.02 

5.75 

4.31 

5.37 

4.33 

5.72 

4.11 

LX  5 

5.91 

3.82 

5.82 

4.00 

LX6 

6.33 

4.31 

6.09 

4.75 

6.90 

4.79 

6.71 

4.65 

6.17 

4.46 

LX7 

5.68 

4.19 

6.69 

4.46 

6.01 

4.24 

LX8 

6.14 

4.28 

5.68 

4.12 

6.59 

4.46 

6.33 

3.93 

5.76 

2.98 

LX9 

6.28 

4.14 

6.85 

4.62 

6.39 

4.19 

/0  ■If' 


Ot  more  than  incidental  interest  are  the  hypotheses 
listed  below. 

In  general  arc  earthquake  and  explosion  mean 
vectors  equal? 


II: 


i  xxi  = 

j-1 


6 

l 

j-l 


Arc  earthquake  and  explosion  means  equal  on  a 
station-by-station  basis? 


(19) 


II: 


1 


6 


(20) 


This  also  gives  simultaneous  confidence  intervals  for 
contrasts  comparing  earthquakes  and  explosions  on  a 
station-by-station  basis. 


Mow  significant  is  the  linear  contrast  M  -  1.9 m 

S  m] 

on  a  s tat ion- by- station  basis? 


II:  -1.9V 


ljl  *  V1 j  2 


1,9V2jl  +  V2 j  2 ' 


(21) 


For  the  input  to  program  BMD-X-63, Dixon  (1970),  we  note 
that  the  12  x  2  parameter  matrix  B  is  given  by 
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_vlll 

V 

v  1 1 2 

V* 
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• 

• 

• 

• 

• 

• 

• 

• 

• 

Vlol 

V 

V162 

= 

l< 

O' 

V211 

V 

v  2 1 2 

-  (V4 
>1 

• 

• 

• 

• 

• 

• 

• 

• 

• 

V261 

V262 

_ 

vf 

V26 

(22) 


The  123  x  12  matrix  X  has  in  each  row  a  vector  which 
generates  an  observation  in  Table  I.  For 
example,  the  component  which  generates  the  observation 
at  PG  for  the  16th  quake,  say  (4.94,  4.90),  should 
have  mean  so  that  the  row  vector  for  X  should  be 
(0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0)  and  similarly  for 
the  other  observations  in  the  table.  In  order  to  generate 
hypotheses  (19)  through  (21),  it  is  sufficient  to  note 
that  they  are  of  the  form  ABC'  =  D  where  (19),  for 
example,  has 
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w> 


A  =  (1,  1,  1,  1,  1,  1,  -1,  -1,  -if  -if  -lt  .i()( 


,  and  P  *  {Ot  Oj , 


The  generalized  linear  hypotheses  program  yielded 
the  estimate  for  B  given  in  Table  II.  The  table  is 
in  the  format  of  equation  (22).  The  table  tends  to 
indicate  that  there  are  some  significant  differences 
between  the  means.  The  standard  errors  are  computed  by 
noting  that  the  program  produces  an  estimate  for  the 
covariance  matrix 


I 


(23) 


and  that  the  variance  covariance  matrix  of  the  rows  of 
is  given  by  £fl(x'x)  *  where  8  denotes  the  Kronecker 
product,  and 


-1 


(X'X)  -  Diag  (.06,  .10,  .07,  .06,  .06,  .09, 

.17,  1.00,  .14,  .11,  .12,  .17). 


(24) 
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TABLE  II 


Estimated  Mean  Mg,  mb  Vectors 
For  Earthquakes  and  Explosions 

Station  m,  M 

Earthquakes  b  s 


NP 

WH 

PG 

RK 

HN 

KN 

Presumed 

Explosions 

NP 

WH 

PG 

RK 

HN 

KN 


5.90(.02)*  4 . 76  C. 04) 

5.00  (.04)  4.83(.06) 

5.05(.03)  4. 74  ( . 04) 

5.05(.02)  4.66(.04) 

5.19  (.02)  4.81(.04) 

5. 05  ( . 04)  5. 10(. 05) 

6. 19  (.07')  4.20  (.10) 

6. 25  (. 40)  3. 76  (. 60) 

5.86(.06)  4.25(.08) 

6. 28 (. 04)  4 . 38 (. 07) 

6.03 (.05)  4.27 (.07) 

5.91(.07)  4 . 35 ( . 10) 


*  One  standard  deviation  of  the  mean. 
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The  table  indicates  that  the  estimates  for  mu  have  a 

b 

smaller  variance  and  that  WH  has  a  large  variance  associated 
with  its  estimated  mean. 

The  test  of  the  hypothesis  that  the  overall  earth¬ 
quake  and  explosion  vectors  are  equal  yields  a  value  for 
the  Wilks  generalized  variance  ratio  U2  1  n2  =  .335 
which  can  be  compared  with  the  II  significance  point 

^2  1  110^'^^  =  *911,  so  that  the  hypothesis  is  rejected 
in  accordance  with  expectations 

The  comparison  of  earthquakes  and  explosions  on  a 
station  by  station  basis  generates  U2  n2  =  .256  as 
compared  with  U2  10Q  (.01)  =  .772.  The  values  of  the 
contrasts  (i.e.  the  difference  between  the  earthquake 
and  explosion  means  on  a  station  by  station  basis)  are 
the  most  useful  and  we  give  simultaneous  951  confidence 
intervals  for  these  contrasts  in  Table  III.  The  95% 
confidence  interval  for  any  contrast  of  the  form  a'Bc 
where  a  is  q  x  1  and  £  is  p  x  1,  is  given  by  (Morrison, 

1970) 


a’Bc  ± 


F  1 

a:q,N-q' 


1/2 


(25) 


where 


Q2  =  (£'  1  c)  (a’ (X’X)_1a) 


(26) 
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TABLE  III 


Estimated  Contrasts  Between  Earthquakes 
and  Presumed  Explosions 


Difference  by  Station 


mb 

M 

s 

M  -1 
s 

*  9mb 

NP 

1 

• 

04 

o 

1  + 

• 

o 

00 

* 

.56 

+  .13 

1.12 

+  .13 

WH 

-.125  +  .19 

1.07 

+  .29 

3.46 

+  .30 

PG 

-.81  +  .07 

.49 

+  .13 

2.03 

+  .13 

RK 

1.22  +  .07 

.28 

+  .11 

3.61 

+  .11 

HN 

1 

• 

00 

04 

1  + 

• 

O 

00 

.54 

+  .12 

2.12 

i  *13 

KN 

-.86  +  .09 

.75 

+  .13 

2.38 

+  .14 

*  951  Confidence  Interval  of  the  mean. 


The  values  for  the  two  terms  in  (22)  can  be  read  from 
the  BML-X-63  printout. 

From  Table  III  we  see  that  all  contrasts  are  sig¬ 
nificantly  different  from  zero.  An  alternate  proposed 
discriminant  is  the  linear  function  represented 
symbolically  as  M  -  1.9  m^.  The  values  of  this  contrast 
are  shown  in  Table  III  along  with  951  confidence 
intervals.  Again  we  see  that  all  contrasts  are  sig¬ 
nificantly  different  from  zero. 

The  fact  that  the  mean  (M  ,  m^)  point  for  the 
earthquakes  is  different  from  that  of  the  explosions 
is  not,  of  course,  solely  indicative  of  discrimination 
capability  since  the  means  of  either  group  can  be 
shifted  anywhere  along  its  Mg  ,  m^  trend  by  suitable 
selection  of  large  or  small  events.  However,  examina¬ 
tion  of  the  Ms,  m^  plots  for  the  data  shows  that  in 
fact  the  separation  in  the  means  is  truly  due  almost 
exclusively  to  the  separation  of  the  trend  lines.  The 
fact  that  the  earthquake  trend  line  tends  to  have  a 
slope  between  1.6  and  1.9  leads  to  the  discriminant 
proposed  above. 

Finally,  we  note  that  the  main  use  of  the  estimates 
for  the  generalized  mean  vectors  (22)  and  covariance 
matrix  (23)  will  be  in  estimating  the  discriminant 
vectors  in  (22)  to  be  applied  to  each  station 
(Ms,  m^j  pair.  The  estimated  discriminant  vector, 

ij  ■  (iij  -  V’  r1  (28) 
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can  be  applied  to  a  data  sample  of  Ng  observed  station 
vectors  to  obtain  the  value  (1) 


d  C  X  2  » •  •  •  X  ^  ) 
s 


=  I  d!  X. 

j  =  l~J 


I (CT)  . 
j  J 


(29) 


to  be  compared  with  a  threshold  C(a)  which  produces  a 
desired  false  alarm  probability  a.  Of  course,  the  rate 
a  is  only  obtained  in  theory  when  the  mean  vectors  and 
covariance  matrices  are  known  exactly.  The  explosion 
detection  probability  will  be  estimated  by  (8)  with 


(30) 


2 

in  place  of  Dj,(5).  Since  the  explosion  detection  prob¬ 
ability  is  an  increasing  function  of  JD?,  ic  is  clear 

j  J 

2 

that  Dj  measures  the  contribution  to  detection  prob¬ 
ability  of  station  j.  The  overall  detection  probability 
can  be  increased  either  by  an  event  which  records  at 
many  stations  or  by  an  event  which  records  at  a 
relatively  few  stations  with  high  D?  values. 

The  above  results  can  be  illustrated  using  the 
data  in  Table  I.  A  discriminant  analysis  program  (DISNP) 
which  accepts  the  generalized  mean  vectors  and  covariance 
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matrix  along  with  each  event  and  the  identities  of  the 

stations  recording  each  event  was  used.  Table  IV  shows 

the  resulting  estimated  discriminant  vectors,  correc- 

2 

tion  terms  and  D  values  for  each  station. 

We  note  that  station  WH  adds  the  most  to  the  detec¬ 
tion  probability  with  NP  adding  the  least.  (Wil  may  be 
good,  however,  mostly  because  of  its  small  data  sample, 
see  Table  I,  an  effect  pointed  out  by  Shumway  and 
Blandford  (1970)).  In  order  to  determine  performance 
measured  by  plotting  detection  probability  as  a  function 
of  false  alarm  rate  (Figure  1)  we  use  equations  (9)  and 
(10)  for  D  =1./,  5.0,  10.0,  and  15.0  in  the  program 
IDCURVE.  All  that  is  necessary  to  evaluate  the  perfor¬ 
mance  of  a  particular  sub-group  of  stations  is  to  add 
up  the  D  values  for  the  sub-group  and  refer  to  the 
identification  curve  in  Figure  2.  We  see  that  the  only 
politically  unacceptable  performance  occurs  when  NP  is 
the  only  station  recording.  Any  two  stations  (exclusive 
of  NP)  which  record  give  a  D2  of  at  least  10  yielding  a 
detection  probability  of  .90  for  a  false  alarm  prob¬ 
ability  .001.  The  worst  combination  of  stations 
recording  (PG,  HN)  in  Table  I  still  yield  a  D2  of 
greater  than  10  or  a  detection  probability  of  greater 
than  .90  in  the  allowable  false  alarm  region.  It  is 
important  to  emphasize  that  the  data  set  used  consists 
of  very  large  events  and  does  not  include  any  "anomalous" 
earthquakes  such  as  those  discussed  by  Der  (1972)  . 

The  values  of  the  discriminant  vector  applied  to 
the  20  earthquakes  and  9  explosions  are  shown  in  Table  V. 

If  we  arbitrarily  assign  earthquakes  equal  prior  odds, 
the  threshold  is  0  and  no  misclassif ication  result. 
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TABLE  IV 

Discriminant  Analysis  Results 


STATION 

discriminant 
“  Vector 

CORRECTION 
lfERM  — 

I)^ 

NP 

(-S.02,  3.43) 

-14.93 

1.70 

WH 

(-16.37,  9.25) 

-52.30 

15.24 

PG 

(-9.81,  5.13) 

-30.48 

5.24 

RK 

(-13.03,  5.84) 

-47.41 

8.79 

HN 

(-10.18,  5.40) 

-37.55 

5.68 

KN 

(-11.25,  6.38) 

-31.41 

7.22 

llJr 


Figure  2.  Identification  curves  for  selected  values  of  I) 
multiple  station  discrimination. 


TABLE  V 


Values  of  Discriminant  Function  (5.10) 


For 

Data 

in  Table 

!  I 

—  —  y 

Events 

NP 

WH 

PG 

RK 

HN 

KN 

Discriminant 

Function 

EQ1 

X 

X 

X 

2.50 

EQ2 

X 

X 

X 

10.98 

EQ3 

X 

X 

X 

7.42 

EQ4 

X 

X 

X 

15.06 

EQ5 

X 

X 

X 

X 

34.15 

EQ6 

X 

X 

X 

X 

17.68 

EQ7 

X 

X 

X 

19.87 

EQ8 

X 

X 

X 

X 

X 

S5.68 

EQ9 

X 

X 

X 

X 

X 

X 

56.76 

EQ10 

X 

X 

X 

21.40 

EQ11 

X 

X 

X 

X 

13.70 

EQ12 

X 

X 

X 

X 

X 

X 

34.45 

EQ13 

X 

X 

X 

X 

10.23 

EQ14 

X 

X 

X 

X 

X 

19.86 

EQ15 

X 

X 

X 

X 

X 

X 

51.69 

EQ16 

X 

X 

X 

X 

X 

X 

55.20 

EQ17 

X 

X 

X 

X 

X 

27.45 

EQ18 

X 

X 

X 

34.97 

EQ19 

X 

X 

X 

X 

X 

X 

63.68 

EQ20 

X 

X 

X 

X 

X 

35.70 

EX1 

X 

X 

X 

-21.03 

EX2 

X 

X 

X 

X 

X 

-18.67 

EX3 

X 

X 

X 

X 

-11.62 

EX4 

X 

X 

X 

X 

X 

X 

-30.83 

EX5 

X 

X 

-12.35 

EX6 

X 

X 

X 

X 

X 

-41.57 

EX7 

X 

X 

X 

-23.17 

EX8 

X 

X 

X 

X 

X 

-36.29 

EX9 

X 

X 

X 

X 

-34.47 

I7cL 


It  is  interesting  to  compare  the  multiple  station 
discriminant  with  one  derived  only  from  event  means 
(the  usual  approach).  Table  VI  shows  the  event  means. 
In  this  case  Ng  ■  1  and  we  obtain  the  discriminant 
vector  (-15.20,  7.31)  with  correction  term  -52.45.  The 
value  of  D  was  7.246,  which  puts  the  theoretical 
detection  rate  at  .75  for  a  false  alarm  rate  .001. 
Figure  2  shows  the  identification  curve  using  the 
sample  means  and  variances  of  the  discriminant  values 
rather  than  the  theoretical  means  and  variances  based 
on  D  .  Table  VII  shows  the  values  for  the  discriminant 
function. 


TABLH  VI 


Data  on  Measured  Body  Wave  Magnitudes  (mb)  (Xj) 

and  Surface  Wave  Magnitudes  (M  )  (X2)  for  9  Presumed 
Explosions  (pop.  2)  and  20  Earthquakes  (Pop.  1) 


Explosions  (E) 


OBSERVATIONS 

POP2  m.  M 

b  s 


6.04 

4.33 

S.  97 

4.39 

5. 85 

4.35 

5.79 

4.14 

5.87 

3.90 

6.51 

4.49 

5.74 

4.22 

5.98 

4.08 

6.07 

4.30 

Earthquakes 

JSL 

OBSERVATIONS 

POP1  mb 

5.60 

4.25 

5.18 

3.93 

6.31 

6.30 

5.36 

4.49 

5.96 

6.39 

5.  26 

4.42 

5.17 

5.10 

4.75 

4.50 

5.35 

5.49 

5.01 

4.48 

5.27 

4.41 

5.27 

4.69 

4.98 

3.66 

5.11 

3,99 

5.06 

4.58 

5.09 

4.90 

5.15 

4.82 

4.56 

4.08 

5.00 

4.94 

5.43 

5.48 

/?-> 


TABLE  VII 


Values  of  Discriminant  Function  for  Data  in  Tablo  VI 


EQ1 

4.61 

EX1 

-7.71 

EQ2 

2.44 

EX2 

-6.21 

EQ3 

2.58 

EX3 

-4.68 

EQ4 

3.79 

EX4 

-5.30 

EQ5 

8.56 

EX5 

-8.27 

EQ6 

4.80 

EX6 

-13.69 

EQ7 

11.14 

EX7 

-3.96 

EQ8 

13.14 

EX8 

-8.36 

EQ9 

11.25 

EX9 

tr 

• 

GO 

1 

EQ10 

9.04 

EQ11 

4.58 

EQ12 

6.62 

EQ13 

3.50 

EQ14 

3.94 

EQ15 

9.01 

EQ16 

10.89 

EQ17 

9.40 

EQ18 

12.96 

EQ19 

12.55 

EQ20 

9.96 

6.  OTHER  POSSIBLE  APPROACHES,  SUMMARY 

Two  other  approaches  were  tried  for  the  multiple 
station  case.  The  first  was  to  use  the  (M  m.)  pairs 

m  xx  stations  as  a  12-element  vector  with  missing 
observations  filled  in  as  sample  means.  This  tended  to 
produce  a  12  x  12  covariance  matrix  which  was  nearly 
singular.  Since  some  columns  became  nearly  linear  combi¬ 
nations  of  the  others  this  tended  to  produce  a  12  x  12 
covariance  matrix  which  was  nearly  singular.  Thus  the 
performance  of  the  discriminant  on  future  data  as  in 
Shumway  and  Blandford  (1970)  could  not  be  evaluated,  and 
doubt  was  cast  on  the  entire  technique.  A  second  approach 
used  a  non-linear  function  composed  of  the  ratio  of 
generalized  variances  under  the  assumption  that  the  suite 
of  station  observations  originated  from  either  a  bomb  or 
an  earthquake.  While  this  procedure  was  appealing  as  the 
generalized  likelihood  solution,  the  test  statistic  was 
not  normally  distributed.  Detection  curves  based  on  a 
normal  approximation  were  inferior  to  those  for  the 
multiple  station  linear  discriminant. 

The  preceding  calculations  indicate  that  using  a 
multiple  station  discriminant  function  when  the  appro¬ 
priate  data  is  available  will  lead  to  an  average  overall 
improvement  in  the  identification  curves.  Consideration 
of  the  results  in  Shumway  and  Blandford,  coupled  with 
the  realization  that  our  new  technique  is  still  basically 
two-dimensional  instead  of  12-dimensional,  suggests  that 
the  improvement  is  real  and  not  due  to  higher  dimensionality 
of  the  new  technique  working  on  a  small  data  base.  However 
this  possibility  is  not  definitely  excluded  in  this  study 
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which  is  aimed  more  at  outlining  the  technique  and  giving 
an  illustrative  example. 

This  improvement  may  be  largely  due  to  the  proper 
accounting  made  of  "station  effects"  which  would  be  the 
same  for  earthquakes  and  explosions  and  which  could  also 
be  accounted  for  by  the  traditional  method  of  "station 
corrections"  applied  before  calculation  of  the  average 
Ms,  m^  values.  The  present  method,  while  preserving  the 
capability  of  discrimination  when  reports  from  some 
stations  are  missing,  can  in  addition  take  advantage 
of  different  Mg  or  m^  patterns  between  earthquakes  and 
explosions,  i.e.  the  possibility  of  different  radiation 
patterns . 

Several  possible  additional  sources  of  improvement 
need  to  be  investigated.  These  improvements  could  be 
possible  using  more  detailed  assumptions  about  the 
covariance  matrices.  For  example,  it  is  generally 
thought  that  the  scatter  for  the  earthquake  distribution 
exceeds  that  of  the  explosion  distribution  although  the 
performance  of  the  classical  discriminant  is  not  degraded 
significantly  by  this  violation  (Shumway  and  Blandford, 
1970).  In  any  case,  a  class  of  linear  functions  of  the 
form 


V  -  <vij  -  V’di  -  'W1 

with  ^2  the  explosion  and  earthquake  covariance 
matrices  and  8  a  parameter  chosen  to  optimize  the 
detection  probability  might  make  a  significant  improve¬ 
ment.  A  second  possibility  is  to  take  advantage  of 
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inhomogeneities  in  the  covariance  matrices  for  different 
stations  by  replacing  l  in  equation  (2)  by  . 

A  problem  of  practical  interest  arises  when  a  few, 
but  not  all,  discriminant  measures  are  missing  from  a 
station.  In  this  case,  one  still  knows  the  marginal 
distribution  of  the  variables  present  to  be  normal  and 
the  analysis  proceeds  as  before.  Let  £  be  partitioned 
into 


11 

hi 

21 

^■ll 

where  all  values  corresponding  to  the  second  set  of  p2 
variables  are  missing.  Then  the  marginal  density  of  the 
first  set  of  P;L  present  variables  is  normal  with 
covariance  matrix  Hence,  in  the  discriminant 
function  (2)  just  replace  l  by  £n  and  0^ .  -  V-.)  by 
(Vlj  "  v2j)  where  the  starred  mean  difference  is^a  px 
by  1  vector  composed  of  the  first  px  components  of 
(vlj  “  v2j) •  In  these  cases,  one  still  gets  the  use  of 
a  station  value  as  a  possible  discriminator  even  though 
it  is  a  vector  of  lower  dimension  than  one  might  prefer 
under  ideal  conditions.  This  requires  that  one  retain  a 
different  covariance  matrix  for  each  station  but 
preserves  the  optimality  properties  that  would  be 
destroyed  by  guessing  values  for  the  missing  variables. 
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